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ABSTRACT 



To probe the weakly non-linear regime, past the point where simple linear theory is 
Q^, sufficient to describe the statistics of the density distribution, we measure the skewness 

{S3) and kurtosis (^4) of the Count Probability Distribution Function (CPDF) of 
the IRAS 1.2 Jy sample obtained from counts in cells. These quantities are free 



method yields unbiased results. 



(SJ ' parameters in a maximum likelihood fit of an Edgeworth expansion convolved with a 

' Poissonian to the observed CPDF. This method, applicable on scales ^ 5h^^ Mpc, 

^ ■ is appreciably less sensitive to the tail of the distribution than are measurements of 

^ ■ 53 and 5*4 from moments of the CPDF. We measure 5*3 and 54 to / ~ 50/i~^Mpc; 

the data are consistent with scale invariance, yielding averages of (^3) = 2.83 ± 0.09, 

CS ! and {S4) = 6.89 it 0.68. These values are higher than those found by Bouchet et 

f~^ I al. (1993) using the moments method on the same data set, {S3) = 1.5 it 0.5 and 

^^ I (S'4) = 4.4 lb 3.7, due to lack of correction in the latter work for finite-volume effects. 

i-Q ' Unlike the moments method, our results are quite robust to the fact that IRAS galaxies 

I ' are under-represented in cluster cores. We use A'^-body simulations to show that our 

cd 
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1. Introduction 

Many approaches have been made to characterize the clustering of galaxies, especially over 
the past decade as better and deeper redshift surveys have become available (cf., Borgani 1995; 
Strauss & Willick 1995, for comprehensive reviews). The two-point correlation function and its 
Fourier Transform, the power spectrum, the long popular methods for describing the clustering of 
galaxies, are complete statistical descriptions of the density field only if the phases of the Fourier 
modes of the density field are random. Indeed, simple inflationary models predict these random 
phases; if this condition holds, the one-point probability distribution function (PDF) of the density 
field 5{r) is Gaussian. 
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As perturbations grow by gravitational instability, an initially Gaussian distribution remains 
Gaussian as long as the fluctuations remain in the linear regime. However, once non-linear efi^ects 
become important, the distribution deviates from its initial Gaussian state, and one needs higher 
order statistics to characterize the density field. 

For a zero-mean Gaussian distribution all reduced moments (cumulants) of the PDF are zero 
except the variance ((5^) = o"^), hence non-zero skewness ((5^), kurtosis ((5^) — Su^, and higher 
order cumulants are measures of the deviation of the distribution from a Gaussian. In this paper, 
we only consider the first two lowest order effects, the skewness and the kurtosis. These A^th-order 
cumulants are equal to the volume averaged correlation functions, 

Cn{v) = nv / d^^id^^2 ■ ■ ■ d^rN CN{ri,r2, . . . ,rN), (1) 

where 1^2 = c"2 = {5"^) , ^^3 = ((5^), ^^4 = {5^) — Su^ and so forth. Here the volume v over which the 
^7V are averaged is defined by the smoothing scale of the density field 5. 

The assumption of scale invariance by Balian & Schaeffer (1988, 1989), 

eiv(Ari, . . . , Ar^) = X'^''-^Hn{vi, • • • , r^), (2) 



yields the scaling relation, 



iN{v)=SNe2- {y)^ (3) 



where the S^ are independent of scale. The scale invariance of Sn and the scaling relation 
@ are in fact predicted by perturbation theory in the mildly non-linear regime, under the two 
assumptions of Gaussian initial conditions and growth of conditions via gravitational instability 
(Fry 1984ab, Bernardeau 1992). Thus one can test the scale- invariance model by measuring 
the dependence of Sn on smoothing scale, although it can be difficult in practice to rule out 
non-Gaussian models (Fry &; Scherrer 1994; Bouchet et al. 1995). 

Calculation of Sn from gravitational instability invokes {N — l)th order perturbation theory. 
Bernardeau (1994b) presents a method for calculating Sn for top-hat filters; the results are: 

53 = f-(n + 3), (4) 

and 

where n is the spectral slope of the power spectrum; a pure power-law spectrum is assumed. The 
expressions above are for Oq = 1; while the Sjy are sensitive to the slope of the power spectrum, the 
dependence on Qq is quite weak (Bouchet et al. 1992; 1995). These, and the corresponding results 
for a Gaussian filter have been confirmed with A-body simulations (Juszkiewicz, Bouchet, & 
Colombi 1993; Bernardeau 1994b; Lokas et al. 1995; Juszkiewicz et al. 1995). On the observational 
side, calculations of ^3 and 5*4 have been done for the CfA (Huchra et al. 1983) and SSRS (da 



Costa et al. 1991) catalogs by Gaztanaga (1992) and Fry & Gaztafiaga (1994), and on the IRAS 1.2 
Jy sample (Fisher et al. 1995) by Bouchet et al. (1993, hereafter B93). Calculations of higher order 
angular moments have been done for the Lick galaxy counts (Szapudi, Szalay, & Boschan 1992), 
IRAS galaxies (Meiksin, Szapudi, & Szalay 1992), the APM galaxy survey (Gaztanaga 1994, 1995; 
Szapudi et al. 1995), and the EDSGC (Szapudi, Meiksin, & Nichol 1996). For optically selected 
galaxies, Gaztafiaga (1992) found that 5*3 = 1.94 it 0.07 up to a smoothing scale of ~ 22 h^^ Mpc, 
which is slightly higher than the value found by B93 for the IRAS sample: 5*3 = 1.5 it 0.5. 

The standard technique to measure 5*3 and 5*4 from observational data involves calculation 
of the moments of the Count Probability Distribution Function (CPDF)^] (§ 2J.), which in turn 



is determined via counts in cells. The high-order moments of the CPDF are of course weighted 
heavily by its high density tail. Regions of such high density are rare, so these moments are highly 
sensitive to the presence or absence of a few clusters (Colombi, Bouchet, & Schaeffer 1994, 1995; 
Szapudi &: Szalay 1996). Also, in a finite volume, there is always a densest region, and thus the 
CPDF goes to zero for higher densities. Not taking this finite-volume effect into account will 
cause the clustering amplitudes to be systematically underestimated. The CPDF asymptotes to 
an exponential at high densities, especially in the strongly non-linear regime (Balian & Schaeffer 
1989); one can thus extrapolate the observed CPDF to arbitrarily high densities. Thus one can 
obtain unbiased estimates of the high-order moments, if the volume is large and dense enough to 
reach this asymptotic regime, to allow the exponential to be fit (Fry &; Gaztanaga 1994; Colombi 
et al. 1994, 1995). This exponential asymptotic behaviour is not expected for the weak regime, 
making it difficult to correct for finite- volume effects (cf., the discussion in B93). Finally, the tail 
of the distribution is also affected by finite- sampling effects; it can be underestimated if the CPDF 
is determined from too few spheres (Szapudi & Colombi 1996). 

However, as we shall see, skewness and kurtosis affect the entire CPDF, not just the tails. In 
particular, skewness causes the mode of the distribution (the region where the measured CPDF 
is most robust) to shift from the mean. This motivates us to develop a new method to measure 
the Sn from the CPDF, less sensitive to finite volume effects, by fitting the entire CPDF to a 
functional form. 

There are several approaches to calculating the evolution of the PDF of 6 from Gaussian 
initial conditions, using the Zel'dovich (1970) approximation (Kofman et al. 1994) or Eulerian 
perturbation theory (Bernardeau 1992; Bernardeau & Kofman 1995; Colombi et al. 1997). The 
so-called Edgeworth expansion, which provides a convenient parametric form to account for small 
deviations from Gaussianity, gives an excellent fit to the PDF of 5 in A^-body models for small 
a (Juszkiewicz et al. 1995). In this paper, we take the Edgeworth expansion convolved with a 
Poisson distribution as our model, and perform a maximum likelihood fit with respect to the 
free parameters Sn to the observed CPDF of the IRAS 1.2 Jy survey from B93. We expect this 



^The PDF refers to the distribution function of the underlying continuous density field S, while the CPDF is the 
distribution function of the discretely sampled galaxy distribution. 



method to be more robust than direct calculation of the moments, since it depends more on the 
overall shape of the distribution function than on the high-density tail region. Although there 
are several applications of the Edgeworth expansion to measure non-Gaussian statistics of the 
density field in the literature (Scherrer & Bertschinger 1991; Amendola 1994; Juszkiewicz et al. 
1995), it has not yet been applied to observational data. We check the validity of our technique 
by applying it to IRAS mock catalogs taken from A^-body simulations, and compare the results 
with the predicted value of 5*3 and 5*4 from perturbation theory, and with values measured from 
the moments of the PDF. 

In §2 we give a brief account of the moments method, and describe our model based on the 
Edgeworth expansion. We test this method with A^-body simulations in §3. In §4, we apply our 
method to the IRAS 1.2 Jy CPDF and compare our results with those of the moments method 
and perturbation theory. We summarize our results in §5. 



2. Method and Analysis 

2.1. Count Probability Distribution Function and its Moments 

The CPDF Pn{1) is defined as the fraction of randomly positioned spheres of radius I 
containing exactly A^ galaxies for a given volume-limited galaxy sample. Here we use the CPDF 
of IRAS galaxies from 10 volume-limited subsamples as calculated by B93 (see their Table 1). We 
place 10^ random points in each subsample and count the number of galaxies in concentric spheres 
of different radii I, from each point, considering only those spheres that are completely included 
within the subsamples (see B93 for details). 

B93 calculate the normalized cumulants Sjy by the moments m,ethod. The moments of the 
distribution P/v(/) are given by 




^m(0= — ^ =>. — ^ Pn{1), (6) 



where A^ = (A^) = J2^Pn{1) is the mean number of galaxies in a sphere of size I. The first few 
volume-averaged correlation functions (reduced moments), corrected for shot noise are given by 
(Peebles 1980) 

6(0 = //2-i, (7) 

6(0 = /^3-3f + i„ (8) 

6(0 = /.4-6|-3^i + ll^-i3. (9) 

The skewness and kurtosis, S3 and S4, then follow from equation (^. These calculations were 
done by B93 for the IRAS redshift survey, and it was found that the scaling relation (equation |3|) 



indeed holds very well (see Figure 8 of B93). A fit of the data to log^A? = Cat + Dn logC2i gives 
D3 = 1.96 lb 0.06 and D^ = 3.03 ± 0.18, where scale invariance predicts Dn = N — 1 (equation |3[). 
All calculations are done in redshift space, but the Sj\f are quite insensitive to redshift space 
distortions, at least on mildly non-linear scales (Bouchet et al. 1992; Lahav et al. 1993; Fry & 
Gaztanaga 1994; Hivon et al. 1995). As mentioned in the previous section, no correction for 
finite-volume or finite-sampling effects have been carried out for these data. 



2.2. The Edgeworth Expansion 

The primordial density fluctuations are assumed to be Gaussian distributed, and as these 
fluctuations grow by gravitational instability, the PDF of 5 deviates away from its initial Gaussian 
form, generating non-zero higher order moments. To the extent that the deviations from a 
Gaussian are small, it makes sense to write the PDF as an expansion around a Gaussian. The 
Edgeworth expansion is a rigorous way to do this, as described by Juszkiewicz et al. (1995). 

We expand the PDF, here denoted by p{i^), where 5 = (p — p)/p, v = (5/(T, in terms of a 
Gaussian 

cl){v) = ^=exp(-i/V2) (10) 
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and its derivatives. This is called the Gram-Charlier series (Cramer 1946), 



p{y) = (j){y) 



co-^^H^{u) + ^H2{v) + 



(11) 



where the Hi are the Hermite polynomials, as given in Table 1. By the orthogonality of the Hi 
one obtains, 

/oo 
Hi{u)p{u)du. (12) 

-00 



Therefore the first few coefficients of equation (|ll|) are given by: 

co = l, ci = C2 = 0, CI = {-l)^ Sia^~^ (3<Z<5), cg = Sgcr^ + Sg V (13) 

where the Si are the normalized cumulants defined in equation (^) . A reordering of the terms of 
the Gram-Charlier series gives a proper asymptotic expansion in a, the Edgeworth series: 

p{u) = (^(i/) |l + j^S^H:,{u)a + [^54^4(^) + ^SIh^{v)^ a^ + 0{a^)^ . (14) 



2.3. A Model for the CPDF 

Equation ( |l^ ) is a model for the underlying density distribution from which the galaxies are 
sampled, but it does not yet account for the discreteness of the galaxies. Due to the finite number 
of galaxies in each sample, the observed CPDF is subject to Poisson noise, and we take this effect 



into account by convolving the Edgeworth expansion with a Poisson distribution (cf., Coles & 
Jones 1991). 

We define the density contrast as 6 = 7^> ' . Let us rewrite the Edgeworth expansion to 
third order (second order in a) as a function of 6, 

E{6) 



7i^-"^^"{^n^-^3(W+[|^^ 



Z4(V^) + ^He{5/a) 



a- . (15) 



The expectation value of P{N) at a given value of / is given by the convolution of equation (15) 
with a Poisson distribution, 



(P(iV)) = j d5 E{5)F{N I Ns = (N) {S + 1)) (16) 

where the Poisson distribution is 

F(Ar I iV^) = ^e-^^ (17) 

the conditional probability of finding N points in a sphere when the true overdensity in that 
sphere is 6 = {Ns - {N))/{N). 

The Edgeworth expansion is not positive definite; moreover, for values of a approaching unity, 
it shows unphysical oscillations (Juszkiewicz et al. 1995; Ueda & Yokoyama 1996). However, 
when we convolve it with a Poisson distribution, the Edgeworth expansion becomes much better 
behaved. Figure |^ illustrates this with the observed CPDF for a volume-limited sample of IRAS 
galaxies (points). The CPDF is normalized to the total number of spheres which fill the volume, 
M (equation 18), and the error bars are given by Poisson statistics, i.e., the square root of the 



value of the CPDF. The solid curve is the best fit of the model in equation (|16| ) (using the method 



described below), while the dashed curve shows the underlying Edgeworth expansion (equation 15) 
with the same values of a, 5*3, and 5*4. The dashed line goes negative and oscillates around the 
CPDF, while the solid line traces it nearly perfectly. This example has a = 0.77. The dotted line 
is the result of performing a fit of the Edgeworth expansion with no Poisson noise term included. 
Not surprisingly, the a is overestimated, and the model gives a very poor fit. Thus for the sparse 
IRAS data we use in this paper, the Poisson noise term is absolutely essential in our model for the 
CPDF. Note that the Edgeworth expansion is no longer valid for cr > 1, and thus the example in 
Figure |l| represents the smallest scales on which we will apply it. 

The observed CPDF, P(N), and the model, {P{N)), are defined to be the probability that a 
sphere of size I contains A^ points. When we assume that a total of M spheres are placed randomly 
within the volume, then the number of spheres that contain exactly A^ points is P{N) = MP{N); 
we similarly define {P{N)) = M{P{N)). We would like M to represent the number of statistically 
independent spheres, in order to allow us to define error bars on the observed P{N), but it is not 
clear a priori how to measure this. The number of spheres thrown within the volume, 10^, is 
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Fig. 1. — Comparison of the fit of the Edgeworth expansion convolved with a Poissonian (solid line) 
to the IRAS CPDF of sample size R = 59/i~^Mpc and smoothing length / = 7.92/i^^Mpc (dots), 
and the corresponding underlying density field (the pure Edgeworth expansion without shot noise 
component; dashed line). The top panel is a linear plot, while the bottom has a logarithmic y-axis 
to better show the tail. The best-fit values of the parameters are a = 0.77, ^3 = 2.56, S4 = 7.38. 
The dotted line is the result of a pure Edgeworth expansion fitted to the CPDF, without convolution 
with a Poissonian; the resulting parameters are a = 0.91, 5*3 = 2.92, 6*4 = 12.95. 



clearly an overestimate for M, due to severe overlap between sphere^. One possibility, which is 
used to define errorbars in Figures |l] and |2|, is to take 

M = ^m'. (18) 



47r V I 

the ratio of the volumes of the sample and the sphereQ; here R is the radius of the subsample, and 
uj is the solid angle it subtends. However, Gaztahaga & Yokoyama (1993) show (as we confirm in 
§ |2.4| ) that this is an underestimate, and suggest multiplying equation ( p!8|) by a~^ (which is much 
larger than unity on large scales). This issue is further discussed in Szapudi & Colombi (1996) 
and Szapudi et al. (1996). We do not have a rigorous solution to this problem. Our approach for 



the present paper is to use the value of M given by equation (18), demonstrate directly that it is 
an underestimate, and, in the following subsection, suggest an empirical rescaling to allow us to 
define error bars on measured quantities. 

As explained above, we throw 10^ spheres in each volume to measure P{N). This is multiplied 
by M, which is several orders of magnitude lower than 10^, to obtain P{N). We assume that 
this number of spheres P{N) is Poissonian distributed around the true value {P{N)), so the 
likelihood of each observed data point P{N) is given by the Poisson distribution. We assume 
that the values of P{N) are statistically independent (We will see in the following subsection that 
these assumptions appear to be violated, but that we can make a heuristic fix to the likelihood). 
We therefore can express the likelihood function of the observed CPDF as the product of these 
quantities over N, 

In practice, the product extends only over those values of A'^ for which P{N) > 1, as the Poisson 
model breaks down beyond that. Nevertheless, we will see that in many cases the best fitting 
curve continues to trace the CPDF reasonably well even when P{N) < 1 (cf. the right-hand 
panels of Figure |2|). We have assumed a Poisson distribution in P{N); if we had assumed a 
Gaussian distribution, our likelihood function could be expressed in terms of x^- Here, the 
quantity corresponding to x^ is then: 

C = -2lnL = J2^[{P{N)) - P{N)ln{P{N)) + ln{P{N)\)\ , (20) 

N 

and we minimize this quantity instead of maximizing L (equation [T9| ). We use Powell's method 
(Press et al. 1992) to minimize C with respect to a, S3 and ^4; indeed, we do fits to first order 
(i.e., Gaussian), second order (including the term proportional to S^cr) and third order in the 



■'lO" may in fact be too small a number of spheres to avoid finite-sampling effects in the tails of the distribution; 
cf., Szapudi et al. (1996) and Szapudi & Colombi (1996). However, we argue below that the Edgeworth method is 
insensitive to such biases in the tail. 

*Note that M, and accordingly P{N), are not integers. 



Edgeworth expansion (equation |T^, although we focus mainly on the third-order fit in this paper. 
To avoid having the code settle into a local mininiuni, we start the search with reasonable values 
of the parameters, i.e., those we calculate by the moments method. In most cases this converges 
to a fair estimate of the parameters, but in some cases, especially on large scales, where the higher 
order correlations are intrinsically weak, or when the volume of the subsample is very large and 
hence the density of the galaxies are small, calculating the moments yields negative quantities of 
5*3 and 5*4. When this happens, we set the initial values of S3 and 5*4 to zero instead, and minimize 
C again. 

Maximum likelihood fits of first, second and third order to the IRAS CPDF's for representative 
values of R and / are given in Table 2; those with asterisks are plotted in Figure g. The left panels 
of the figure give P{N) on a linear scale to show the fit around the peak clearly. The right panels 
are log plots to show the tail region in detail. The error bars are again given by Poisson statistics 
as in Figure [l|. The solid line is the first order (Gaussian) fit, and the dotted line and the short 
dashed line represent second and third order fits respectively. The second and third order fits are 
a remarkable improvement over the initial set of parameters obtained by the moments method, 
represented by the long dashed lines. The third order fits, with the term proportional to o"^ 
included, are an improvement over the second order fits. Note in particular that the third-order 
fits follow the observed CPDF even into the region where P{N) < 1, which is not even included 
in the fitting procedure. On the other hand, quite interestingly, the third order model with the 
parameters given by the moments method (long dashes) usually does the best job of fitting the 
high-density tail of the CPDF, although this fit is poor at intermediate to low 5. This again 
reflects the sensitivity of the moments method to the tails of the distribution. As the volume 
of the subsample increases, going downward in the figure, the Gaussian fit becomes a better 
approximation to the observed CPDF. Indeed, for the bottom-most panel with a subsample of 
volume 176 h~^ Mpcand a smoothing length of 31.54 h~^ Mpc, the first, second, and third order 
fits are barely distinguishable. This is a result of two effects. First, the larger the subsample is, 
the sparser the sample becomes, and therefore the Poisson noise dominates the CPDF, swamping 
the higher order correlations. Second, the larger the smoothing scale, the weaker the higher order 
correlations become, because clustering is weaker at larger scales. Therefore in either of these 
cases, the accuracies of the higher order correlation terms, ^3 and especially ^4, drop appreciably, 
and the error bars we derive for these quantities in the next section are therefore quite large. 



2.4. Error Estimates and the Value of M 

We determine errors on the parameters from the inverse of the Hessian matrix of the likelihood 
function. Close to the minimum, we expect the likelihood function to be well approximated by a 
quadratic form, 

£(x) ^ £(xo) + -(x - xo) • D • (x - xo), (21) 
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Fig. 2. — Various fits of tlie Edgeworth expansion convolved with a Poissonian to tlie observed 
CPDF (denoted by dots), with error bars given by Poisson statistics. R is the size of the IRAS 
sample, and I the smoothing scale. The long dashed line represents the model using the a, S3, and 
5*4 of the moments method, the solid line our first order fit, the dotted line the second and short 
dashed the third order fit. The volume of the subsample increases as one goes downwards, and the 
solid line (Gaussian) becomes a better approximation. The right panels use a logarithmic y-axis 
to show the tails. The fits trace the data points reasonably well even when P{N) < 1, which were 
excluded in the likelihood analysis. 
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where x is the vector of parameters, and xq is the value of this vector at the minimum value of C. 
Here D, the Hessian matrix, is the second derivative matrix of £. Since the form of the likelihood 
function is given, the Hessian matrix is known to us. The covariance matrix is then obtained by 

[C] = 2D-^ (22) 

and \fCii are the lo" confidence intervals of the parameters Xj. 

The values of the errors derived in this way are of course critically dependent on the value of 
M assumed. The third-order model fits the points in Figures || and g much too well, given the 
Poisson errors shown, suggesting that our value of M is underestimated. We can quantify this 
with a x^ statistic: 

where Poisson error bars are assumed, and the sum is over only the Af values of N for which 
P{N) > 4 (such that the correspondence between the Poisson distribution, and the Gaussian 
distribution assumed in equation (|2^) , is valid). This quantity is less than the number of degrees 
of freedom AA — 4 by factors of as much as ten. As we discussed in the previous subsection, it is not 
clear how to set M a priori, and consequently the errors on P{N) and parameters derived from 
it (cf., the exhaustive discussion in Szapudi & Colombi 1996 for the strongly non-linear regime; 
there does not yet exist a complete theory for the errors in the weakly nonlinear regime). We 
find that the Gaztahaga & Yokoyama (1993) suggestion of multiplying M by a~'^ does not bring 
X^ /{M — 4) to unity, although it does go in the right direction. 

Our approach is an empirical one: we scale our errors of P{N) to force x^ = AA — 4; 
equivalently, we multiply the errors in o", ^3, and 5*4 derived above by {x^ /{M — 4))^'^. This is 
done o posteriori, and so does not affect the best-fit values of the parameters, but it does of course 
have a strong effect on the derived errors. This is not a statistically rigorous procedure, but we 
will justify these errors empirically in the next section. The 1 a errors determined in this way are 
included in Table 2. 

Table 2 includes the values of £ for each fit; these have been scaled by the factor \^ /{M — 4). 
This way one can quantify the goodness of fit of the curves seen in Figure |^ by comparing the 
relative likelihoods. The lower the value of C, the better the fit becomes for given values of R and 
/; a difference of unity is significant at the 1 a level. These values of C are not comparable across 
subsamples or smoothing scales, since in each case the input CPDF is different. 



3. Tests vi^ith A'^-Body Simulations 

To check the validity of our method, we generate IRAS mock catalogs from the A^-body 
simulations of Protogeros & Weinberg (1997), kindly given to us by D. Weinberg. These assume 
an initial power spectrum P{k) oc k~^, and Qq = 0.3, and are evolved forward to the point when 
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as, the rms fluctuation amplitude within an 8 h^^ Mpc sphere, is equal to 0.8. The simulations 
use a staggered mesh PM code by Park (1990), and are evolved within a cube of size 300 h~^ Mpc, 
containing 150^ particles and a 300^ mesh. We assume that the galaxy distribution is unbiased 
relative to the dark matter. We choose a random point within the simulation to represent the 
Local Group, and produce 10 concentric volume-limited mock catalogs with exactly the same 
volumes and number densities as in the real universe (cf., Table 1 of B93). We compute the CPDF 
for these samples with tophat filters of the same radii as are used in the real universe, and fit 
these to our Edgeworth model exactly as was done above. We wish to compare our results to the 
predictions of perturbation theory, and therefore work in real space, not in redshift space. 

For Gaussian filters, A^-body simulation checks for S^ and 5*4 (Juszkiewicz et ai. 1993; 
Bernardeau 1994b; Lokas et al. 1995) have successfully matched the predicted values of S^ and 
S4 for the n = —1 power spectrum. The counts in cells method uses a top hat filter, for which 
perturbation theory predicts S^ = 2.86 and Sf = 13.89 for the n = — 1 power spectrum of our 
simulations, from equation (Q) and (||)^. Unlike the real IRAS data, biasing is not an issue we 
have to consider in the A^-body simulations, and therefore we expect the perturbation theory 
prediction to agree quantitatively with the results of our method. 

Figure H shows the values obtained for S3 as a function of scale from the A^-body model. 
Open squares denote the weighted average of the values obtained by the third-order Edgeworth 
model to the CPDF's of the mock catalogs, at each scale. The data appear consistent with 
scale-invariance, as indeed one expects for this power spectrum (cf., the discussion in Colombi, 
Bouchet, &: Hernquist 1996). Averaging over all scales gives {S3) = 2.93 it 0.09 (dotted line), in 
beautiful agreement with the perturbation theory value (solid line). 

The open circles indicate the results of applying the moments method to the A^-body 
simulation sampled at a density of 0.01 Mpc~. At this sampling density, the Poisson noise is 
relatively small, and the moments method yields an equal-weight average {S3) = 2.90 ib 0.64, 
remarkably similar to the perturbation theory results, and very much consistent with the 
Edgeworth approximation results of the mock catalogs. Note that the scatter around the mean of 
the determinations of 6*3 using the moments method from the densely sampled data is appreciably 
larger than that of the Edgeworth expansion results based on the much sparser mock catalogs. 

When we apply the moments method to the more sparsely sampled mock catalogs (triangles 
in Figure y), we find that the values of 53 are consistently biased low, an effect which worsens 
at larger scales. This is due to the finite-volume effect discussed in § |l], and perhaps also 
finite-sampling effects as well. We have found that the CPDF for these sparse samples on large 
scales never clearly reaches the asymptotic exponential tail discussed, e.g., in Colombi et ai. 
(1995), and thus we cannot fit the tail to correct for these effects. 



^Actually, the values of S^ and S4 we quote are for f2o = 1- For Qq ~ 0.3, Sg — 2.89 (Bouchet et al. 1995), 
differing by only 0.03 from the value quoted above. We are not aware of a direct analytic calculation of the fio 
dependence of 5*4, but Bernardeau (1994a) shows that it is very insensitive to fio- 
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Fig. 3. — Weighted average of 6*3 vs scale, of the Edgeworth approximation apphed to IRAS mock 
catalogs (open squares). The average value {S3) = 2.93 ± 0.09 (la) is indicated by the dotted 
line. The open circles are the results of the moments method applied to an A^-body simulation 
sampled at 0.01 Mpc^^, yielding an average (^3) = 2.90 it 0.64 (Icr). The solid line on the left 
is the perturbation theory prediction S^ = 2.86 for n = —1. All three methods agree well. The 
triangles are the results for the moments method applied to the sparsely sampled mock catalogs. 
It breaks down at scales larger than / ~ 16 /i^^ Mpc due to finite- volume effects (which we do not 
correct for here), and results in a biased value of S3. 
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Fig. 4. — ^4 vs. scale for the mock catalogs. Open squares denote the weighted average of 
the Edgeworth approximation applied to mock catalogs, as in Fig ^. The average over scales is 
{S4) = 13.54 lb 0.53 (dotted line), in agreement with the perturbation theory prediction 5*4 = 13.89 
(solid line) and that of the moments method applied to the full A^-body simulation (open circles). 
The moments method on sparsely sampled mock catalogs (triangles) breaks down on scales larger 
than / ~ 10 h'^ Mpc. 
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Figure Q shows results for ^4; the symbols have the same meaning as in Figure |3|. The 
average value of S4 obtained from the mock catalogs via the Edgeworth approximation is 
(S'4) = 13.54 lb 0.53 (open squares and dotted line), consistent with the perturbation theory 
prediction of S'4 = 13.89 (solid line). The moments method applied to the densely sampled 
A^-body points (open circles) gives (S'4) = 12.82 ± 6.62, which is very similar to the Edgeworth 
approximation results (open squares), albeit again with larger errors, whereas the moments 
method results from the mock catalog (triangles) begin to break down at even very small scales. 
We conclude that our method gives unbiased estimates of S3 and S4, while the moments method 
is systematically biased low for sparse samples by finite-volume effects, which cannot easily be 
corrected for. 

The last test we perform on A^-body simulations is to check the estimation of our error bars, 
via the method explained in § ^^. For several specific values of sample size R and smoothing 
scale /, we draw a series of 50 mock catalogs randomly from the simulation, compute P{N) for 
each, and fit them to find S3 and S4, and their estimated errors, in each case. We found in every 
case that the mean of the estimated errors agreed with the standard deviation of the individual 
values of ^3 and S'4 to within 10%, implying that the error estimates are correct in the mean. Of 
course, this only works when we scale C by the ratio of x^ to the number of degrees of freedom, as 



discussed in § 2.4; if we do not do this, our errors are overestimated by factors of 2 or 3. However, 



we emphasize here that our error estimation is not done in a rigorous way; our approach in § |2.4| is 
empirical at best, and further tests of our errors with simulations with a variety of power spectra 
are needed to justify these errors fully. 



4. Results for IRAS Galaxies 

Table 2 shows that the best-fit values of a, S3 and S'4 are often considerably different from 
those found by the moments method. Figure ^ shows logo" = 1/2 log ^2 versus the smoothing scale 
/, for the values found in the third order Edgeworth expansion fit. The upper panel shows all 
cases separately, with different symbols representing different sample sizes and lines connecting 
the points found within each subsample. As the size R of the subsample increases, a{l) grows for 
any given I; more luminous galaxies show stronger clustering than do less luminous galaxies (B93). 
This luminosity effect agrees well with that of the moments method quantified in B93. However, 
within each subsample the values trace a power law reasonably well, as does the weighted average 
given by the squares in the lower panel. The dotted line in the lower panel is a least square fit to a 
power law. The average of the slopes at each value of R is 7/2 = 0.87 it 0.08, crossing unity at an 
average value of Iq = 5.07 it 1.45/i~^Mpc, while the moments method (triangles in the lower panel) 
yields 7/2 = 0.80 ± 0.03, /q = 5.44 it 0.53/i~^Mpc (B93). This is quite reassuring, considering the 
sensitivity of a to the order of the fit (Table 2). The large error bar in /q is due to its sensitivity 
on luminosity. 

The weighted average values of S'3 and S'4 found by fitting the third order Edgeworth 
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Fig. 5. — A log-log plot of determinations of a(l) vs / from the third order Edgeworth expansion fit. 
In the upper panel, points from each subsample are connected; different symbols denote different 
sample sizes as shown in the legend. Each symbol is used for several samples: the solid line connects 
the points for the smallest sample, the dotted line the next, and the dashed line the largest sample 
that the symbol denotes. The points are slightly staggered to show the error bars. The lower panel 
shows squares as weighted averages of the top panel. The dotted line is a least square fit to the 
squares, and the triangles in the lower panel are an equal weighted average of the results from the 
moments method. 
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expansion are summarized in Tables 3 (averaging at a given sample size R) and 4 (averaging at 
a given scale /), and are plotted in Figures ^, ^, and §. Note that our likelihood method gives 
meaningful error bars on the values of S^ and S4, enabling us to perform a weighted average 
of our determinations on different scales and from different subsamples, assuming that they are 
independent of scale^. 

The figures also show results for the moments method. They are systematically lower than 
the Edgeworth results, due to finite-volume effects. As with the A^-body tests described in the 
previous subsection, we found that on these large scales, the CPDF never reached the asymptotic 
form that would allow us to correct for these effects (compare with Fry & Gaztahaga 1994, 
who did perform this correction, but only on scales below 10/i~^Mpc). Error bars for values 
determined from the moments method can be calculated (Szapudi & Colombi 1996), but are quite 
complicated. Following B93, we use the standard deviation of determinations of 5*3 and S4 as our 
error bar for the moments method, including only positive values. 

The Edgeworth approximation yields values of S3 and ^4 (open squares in Figures ^, |^ and 
^) which show no statistically significant dependence on R or I, consistent with the scale-invariant 
hypothesis (cf., discussion below). Without finite- volume corrections, this scale- invariance is not 
apparent from the moments results. This is not an issue with the Edgeworth method; it is quite 
insensitive to the tail of the CPDF, and thus no correction for finite-volume effects is needed. 

We find a global average of {S3) = 2.83 ± O.OqJ, which is appreciably larger than the value 
{S3) = 1.5 lb 0.5 found in B93. This is consistent with the bias in the moments method we found 
in § ^; we showed there that the Edgeworth method is unbiased, and we therefore believe that 
our determination of (53) for IRAS galaxies supersedes that of B93. Our results are in better 
agreement with the moments methods results for IRAS from Fry &: Gaztahaga (1994), which 
found 53 = 2.1 ± 0.3, ^4 = 7.5 ± 2.1, in the range 3-10 /I'^Mpc, and Meiksin et al. (1992) from the 
angular distribution of IRAS galaxies (53 = 2.2 ± 0.2, S'4 = 10 ± 3 in the range 4-10 h'^Mpc]^. 

Our result uses the third-order Edgeworth expansion; if we use the second-order expansion, we 
find (5*3) = 2.80 lb 0.11, which is consistent with the value from the third order fit. In both cases, 
the average is dominated by the point with the smallest errors, (5*3) = 2.87 ib 0.06 at 8/i~^Mpc. 

Given the estimated errors on 5*3 for each value of R and /, we can do an a posteriori test of 
the hypothesis that 53 is independent of sample, by computing the x^-like statistic: 

n= Y. ^^'''~,^^'^'^\ (24) 

realizations i S3,i 



^There is some covariance due to the fact that the same galaxies are used for the determination of the parameters 
on different scales; we ignore this effect here. 

^ This error bar is smaller than our final quoted value of 0.09 ; see discussion below. 

Note, however, that one does not expect the angular and spatial Sn to be identical (cf., Bernardeau 1995). 
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Fig. 6. — Weighted averages of S3 results from the Edgeworth approximation as a function of 
subsample size R (open squares), from Table 3. The triangles are the results of the moments method. 
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the global average of the Edgeworth approximation results, {S3) = 2.83 it 0.09. 
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where the sum is over all different values of R and /. We find that TZ exceeds the number of 
measurements of S^ by a factor of two. This is due to several effects: cosmic variance, the 
approximation we have made that the errors in S3 determined at different scales from the same 
data are statistically independent, and possible real scale-dependence of ^3 caused by higher-order 
effects (cf., Colombi et al. 1996). Our final error bar should reflect this increased scatter, and so 
we multiply the error on {S3) by TZ^'^ = 1.4, yielding (S^) = 2.83 ± 0.09. It would be interesting 
to separate out these effects with a n = — 2 simulation, in which these higher-order effects should 
be more important (Colombi et al. 1996); this is work for the future. 

Figure |^ shows the determination of 5*4 as a function of scale. For scales larger than 
log/ = 1.3, S'4 becomes very uncertain, with large error bars. However, the weighted average 
values of S4 in the bottom panel stay close to the total average value of (S'4) = 6.89 ± 0.48, with 
error bars overlapping that value at all scales. We carried out the x^ test of equation (^^ for 
S4 as well, finding IZ = 2 again. We thus also multiply our quoted error on (S4) by 1.4, yielding 
(S4) = 6.98 lb 0.68. This again is larger than, albeit statistically consistent with, the moments 
method results of (S4) = 4.4 it 3.7 (B93), and is in agreement with Meiksin et al. (1992) and Fry 
& Gaztahaga (1994), quoted above. 

Gaztahaga (1992) used the moments method to find (S3) = 1.94 it 0.07 for optically selected 
galaxies. Gaztaiiaga (1992) and B93 claim that the lower value of S3 for IRAS galaxies can be 
attributed to the underrepresentation of IRAS galaxies in dense cluster cores (cf., Strauss et al. 
1992). As a test of this, B93 gave extra weight to the IRAS clusters to match the overdensities seen 
in optically selected galaxies; they found the Sat to be quite sensitive to this: (S3) = 3.71 ± 0.95 
and (S4) = 23.6 it 12.1. This demonstrates the high sensitivity of the moments method to dense 
clusters; the moments are heavily weighted by the tail of the CPDF, and this effect is even more 
important for S4 than for S3. This sensitivity is very dangerous, given the fact that the tail is 
generally hard to measure with accuracy, as we have seen. We have fit the Edgeworth model to 
the CPDF of these cluster-boosted counts, and found (S|) = 2.65 ±0.09, (S|) = 7.79 ±0.67, within 
2 standard deviations of the unboosted results above. This is as expected; boosting the clusters 
only affects the CPDF in the tail, and therefore this has only a small effect on our fits. 

One might argue that this result is misleading; if our Edgeworth model is a good fit to the full 
CPDF before the cluster-boosting, it cannot be a good fit afterwards, because the tail has changed 
dramatically, even though it had very little effect on the rest of the CPDF. Clearly, the moments 
method and our method cannot agree in both cases, independent of issues of finite volume effects. 
With our method, S3 and S4 are determined from a fit to that part of the CPDF that is close to 
mean density, and therefore is not highly non-linear, while the moments method is quite sensitive 
to the non-linear tail. Thus our method can measure S3 and S4 in the weakly nonlinear regime, 
even when strong clustering is present. 

The effective power- law index for IRAS galaxies is n = —1.4 (Fisher et al. 1993), which 
would predict that S3 = 3.26 and S4 = 22.4. Why do our results differ from these values? If the 
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IRAS galaxy distribution were biased with respect to the underlying mass, one would expect that 
the skewness and kurtosis would be systematically affected. The linear bias model, 6g = b5M, 
where Sg is the observed galaxy density field and 6m is the underlying mass density contrast, 
predicts S^^g = Ss^u/b and S^^g = S^^M/b'^- However, as Fry Sz Gaztahaga (1993) point out, we 
cannot consider high-order correlations without also considering the possibility of high-order bias 
(Gaztafiaga & Frieman 1994; Juszkiewicz et al. 1995): 

OO r 

^9 = m = E i^M- (25) 

This leads to: 

Ssg = S^M/bi + 2,b2/bl + 0{{5')) (26) 

S^g = S^M/b'' + l2S^Mb2lb\ + Ab^/bl + l2bl/b\ + 0{{5'')). (27) 

Thus without external information on the detailed form of the biasing relation, we cannot make a 
direct comparison of our results with those from perturbation theory. 



5. Conclusions 

We have measured the count probability distribution function via the counts-in-cells method 
for 10 volume limited subsamples of the IRAS 1.2 Jy redshift survey, exactly as in B93. There 
are various approaches to measure the skewness 5*3 and kurtosis 5*4 of the probability distribution 
function of the underlying density field. B93 calculated these quantities for this sample using the 
moments method. They found scale invariance of 5*3 to Z ~ 25 h'^ Mpc, with an average value 
of 5*3 = 1.5 lb 0.5. However, the moments method is very sensitive to the high density tail of the 
CPDF, making the values of 5*3 and S^ sensitive to finite-volume and finite-sampling effects. These 
effects can be corrected for (Fry & Gaztahaga 1994; Colombi et al. 1995; Szapudi & Colombi 1996; 
Szapudi et al. 1996), at least in the strongly non- linear regime with dense sampling. We work 
here with a very sparse redshift survey in the weakly non-linear region (^2 < l), and find that the 
CPDF never properly reaches the asymptotic limit that allows one to correct for finite-volume 
effects. We propose a method less sensitive to the tails of the CPDF, a maximum likelihood fit of 
the Edgeworth expansion convolved with a Poissonian, to the observed CPDF. The Edgeworth 
expansion is valid only in the weakly nonlinear regime (a < 1); unlike the moments method, it 
cannot be applied on very small scales. 

We have tested our method with IRAS mock catalogs extracted from iV-body simulations; 
we find that the derived values of S3 and ^4 are consistent with the analytic predictions from 
perturbation theory, as well as from the moments method as derived from densely sampled 
A'^-body points. Moreover, our estimated errors are consistent with the scatter in 6*3 and 5*4 seen 
in multiple realizations of the sample. The results from the moments method in these sparse 
mock catalogs are systematically biased low, especially on large scales, due to finite- volume effects. 
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Hence we conclude that the Edgeworth approximation is much more rehable and robust than is 
the moments method, especially in sparse samples and in the weakly non-linear regime, where 
there is no simple method to correct for these effects. 

The resulting values of S3 and S4 are found to be {S3) = 2.83 ± 0.09 and (S'4) = 6.89 ± 0.68, 
significantly higher than the results of B93, but consistent with Meiksin et al. (1992) and Fry 
& Gaztahaga (1994). These results are quite insensitive to the fact that IRAS galaxies are 
underrepresented in cluster cores. Both 5*3 and S'4 are independent of scale within the errors from 
5 h-^Mpc to 50 h-^Mpc. 

We have shown that the data are consistent with the scale-invariant hypothesis. It would be 
very interesting to compare these results with those from various specific models with non-Gaussian 
initial conditions, to see at what level we might be able to rule them out. 

Application of the Edgeworth approximation to optical samples should be interesting, 
especially since previous work has shown discrepancies with the IRAS sample, attributed to the 
underrepresentation of IRAS galaxies in clusters. Also interesting would be to apply this technique 
to angular surveys such as the APM galaxy sample, where we could carry this technique out to 
appreciably higher order. We also look forward to applying this technique on the spectroscopic 
and photometric data of the Sloan Digital Sky Survey (cf., Gunn & Weinberg 1995), which will 
allow us to probe appreciably larger spatial scales. 

There is also further work to be done on the method itself. Our understanding of the errors 
and covariances in the CPDF, and therefore the errors in our derived parameters, is poor, and 
thus our final errors are not rigorously justified. As we have seen in § ^, the x^ test (equation 24) 



suggests that our method of obtaining error bars can hide interesting higher order effects. In 
particular, without good a priori errors, we cannot do a proper test of goodness of fit of our 
model. Further analytical work in this direction is needed, together with more extensive tests with 
simulations over a wider range of conditions. 
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of an Alfred P. Sloan Foundation Fellowship, as well as the support of NASA Astrophysical Theory 
Grant NAG5-2882. RSK acknowledges the support of an Assistant in Research Tuition Award 
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Table 1. Hermite Polynomials 





^oH 


= 1 




Hi{v) 


= u 




H2{u) 


= u^-l 




Hsiu) 


= iy^-3u 




H^{u) 


= j,4 _ gj,2 ^ 3 




H^{v) 


= u^ - 10z^3 + 15u 




Heiiy) 


= j,6 _ i^^4 ^ 45j,2 _ 


-15 



Table 2. Results of representative fits 
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Max. Radius^ l^ iV = order a{l) Ss^l) S^il) C "^ 

*37.46 7.92 6.57 m^ 0.70 f.97 3.22 384.49 



f 0.6f±0.09 ... ... 418.27 

2 0.56 ±0.04 4.15 ±0.83 ... 369.02 

3 0.56 ±0.05 3.17 ±0.66 10.82 ± 4.08 363.01 



*46.01 7.92 3.54 m 0.73 1.74 2.69 895.01 



1 0.68 ±0.08 ... ... 1050.8 

2 0.68 ±0.04 2.70 ±0.44 ... 889.81 

3 0.78 ±0.10 2.00 ±0.28 6.80 ± 0.93 882.27 



57.34 7.92 2.28 m 0.85 2.13 6.49 1187.2 



1 0.74 ±0.06 ... ... 1846.9 

2 0.73 ±0.05 3.43 ±0.56 ... 1187.1 

3 0.77 ±0.02 2.56 ±0.11 7.38 ± 0.76 1170.6 



*57.34 12.56 8.39 m 0.53 2.26 4.90 2066.3 



1 0.48 ±0.08 ... ... 2163.1 

2 0.46 ±0.02 3.53 ±0.57 ... 2034.0 

3 0.45 ±0.02 3.15 ±0.32 7.83 ± 3.46 2033.9 



72.35 7.92 0.98 m 0.81 1.30 1.12 701.67 



1 0.67 ±0.06 ... ... 1005.0 

2 0.74 ±0.06 3.64 ±0.79 ... 527.17 

3 0.73 ±0.02 2.85 ±0.11 2.51 ± 1.94 499.30 



72.35 12.56 3.82 m 0.54 1.32 4.35 140.29 



1 0.57 ±0.08 ... ... 142.70 

2 0.57 ±0.02 1.39 ±0.46 ... 138.81 

3 0.52 ±0.12 1.28 ±0.80 -2.4 ± 14.5 138.63 



90.02 7.92 0.45 m 0.87 1.63 1.34 1493.5 



1 0.58 ±0.07 ... ... 1635.6 

2 0.74 ±0.06 4.69 ±0.88 ... 935.40 

3 0.76 ±0.11 3.45 ±0.66 9.17 ±1.70 929.90 



90.02 12.56 1.74 m 0.57 0.86 -3.6 253.75 



1 0.59 ±0.08 ... ... 268.85 

2 0.59 ±0.02 2.01 ±0.56 ... 250.76 

3 0.58 ±0.10 1.61 ±0.73 1.44 ±10.4 250.87 



113.02 19.71 3.17 m 0.37 0.77 -21.6 357.50 



1 0.37 ±0.07 ... ... 353.54 

2 0.37 ±0.03 3.40 ±2.49 ... 348.66 

3 0.33 ±0.06 3.00 ±1.96 -41.36 ± 100.8 347.29 
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Table 2 — Continued 



Max. Radius'^ 


/b 


Arc 


order 


ail) 


^3(0 


s,{i) 


C^' 




*175.99 


31.54 


1.78 


m 
1 
2 
3 


0.34 


-0.38 


-39.66 


551.92 




0.35 ±0.12 
0.35 ± 0.03 
0.40 ±0.16 


1.77 ±0.94 
1.93 ±1.40 


16.29 ±13.6 


546.20 
545.30 

544.77 



* Examples shown in Figure 2 
'^Radius of subsamples (/i^^Mpc) 

Radius of the sphere (/i~^Mpc) 
'^Average number of points in the sphere of size I 

C = — 21nL, L is hkehhood; values shown here are scaled to the x^ P^r degree of freedom (see 



*^Initial values of the fit, obtained by moments method 
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Table 3. S3 and S4 for different subsamples 



Max. Radius (/i^^Mpc) ^3 ^4 



Moments'^ Edgeworth Moments Edgeworth 



23.14 1.15 2.30 ±0.37 0.12 7.09 ±1.25 

37.47 2.11 ±0.14 3.18 ±0.46 10.6 ± 3.87 10.8 ± 9.62 

46.01 1.74 ±0.01 2.31 ±0.24 1.69 ± 1.00 6.59 ± 0.93 

57.34 2.19 ±0.06 2.63 ±0.11 5.70 ± 0.80 7.39 ± 0.75 

72.35 1.26 ±0.06 2.85 ±0.10 0.86 ± 2.96 2.35 ±1.92 

90.02 0.38 ±0.93 3.15 ±0.13 -14.7 ±17.2 8.10 ±1.66 
113.0 0.64 ±1.00 1.98 ±0.67 -12.2 ±25.1 -16.9 ±10.6 
140.3 -2.55 ±1.57 2.60 ± 2.36 -30.5 ± 6.75 -7.61 ± 14.1 
176.0 -8.96 ±10.8 1.73 ± 1.23 -54.2 ± 58.3 17.4 ±11.8 
219.7 -8.71 ±6.91 1.53 ±1.53 -167. ± 178. 14.5 ±43.1 



'^Moments method : Equal weighted average. Errors la dispersion 
Results of Edgeworth expansion fit. Weighted average 
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Table 4. S3 and 6*4 as a function of scale 



I (^-^Mpc) 


^3 




S. 


1 




Moments'^ 


Edgeworth*^ 


Moments 


Edgeworth 


5.00 


1.15 


2.30 ±0.37 


0.12 


7.09 ±1.25 


7.92 


1.64 ±0.37 


2.87 ±0.06 


4.36 ±3.57 


7.06 ± 0.56 


12.6 


0.21 ±2.78 


2.67 ±0.21 


3.18 ±12.7 


5.29 ±1.53 


19.9 


-0.39 ±1.15 


2.30 ±1.19 


-30.5 ±21.8 


14.6 ±12.9 


31.5 


-2.35 ±3.11 


1.02 ±0.76 


-40.2 ±5.92 


-10.0 ±18.0 


50.0 


-16.6 ±9.02 


1.17 ±4.32 


-198. ±161. 


-11.0 ±122. 


All 


-1.76 ±6.27= 
1.54 ±0.47= 


2.83 ± 0.06^^ 


-31.9 ±79.0= 
5.8 ± 7.09= 


6.89 ± 0.48^^ 



'^Moments method : Equal weighted average. Errors la dispersion 

'^Results of Edgeworth expansion fit. Weighted average 

=Equal weighted average of all values 

Weighted average of all values. Note that our final quotes are (53) 
{S4) = 6.89 ± 0.68 (see discussion in § ||). 

=Equal weighted average of all positive values 



2.83 ± 0.09 and 



